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A STUDY OF DAMPING IN NONLINEAR OSCILLATIONS 

By 

Maurice L. Rasmussen 
Stanford University 

and 

Donn B. Kirk 
Ames Research Center 


SUMMARY 


An investigation is made of nonlinear oscillations in which the damp ing 
and static moments are represented by arbitrary polynomial functions of the 
dependent variable. When the nonlinear damp ing is small but the static nonlin- 
earities arbitrarily large, an approximate solution is established which leads 
to expressions for the damping decrement involving elliptic integrals and 
gamma, functions in special cases. An "effective linear damping" is defined and 
a generalized formula for this parameter is obtained that is valid for a wide 
range of nonlinearities in both the damping and static moments. This formula 
is useful, for instance, in deducing the dynamic-stability parameters of 
missiles observed in nearly planar motion in free flight. 


INTRODUCTION 


Many physical systems are describable only in terms of nonlinear ordinary 
differential equations . One example is the large -amplitude pitching motion of 
hypersonic aerodynamic configurations for which the frequency can be strongly 
dependent on amplitude . 1 In the analysis of ballistic -range data, it is 
desired to deduce the aerodynamic properties, linear and nonlinear, from obser- 
vations of a given set of oscillations. Toward this end, a knowledge of the 
effects of the nonlinearities involved in the differential equations governing 
the motion is very important. 

In this study, we will be concerned chiefly with the damping of oscilla- 
tions of nonlinear systems. Although we will be interested in solutions to the 
nonlinear equations for their own sake, the eventual goal will be to derive 
formulas allowing the extraction of both linear and nonlinear damping param- 
eters from a set of data. From the standpoint of aerodynamics, this goal is 
equivalent to determining from data the dynamic- stability parameters of a mis- 
sile oscillating in a plane . A theoretical study indicating the form of the 

••Although the present analysis is applicable to a wide variety of oscil- 
latory phenomena, the text is cast in terms of pitching oscillations of bodies 
having nonlinear damping and static moments in flight in a medium of constant 
density. 


nonlinearities to be expected in the dynamic- stability parameters can be found 
in reference 1. Aspects of obtaining free-flight damping in a conventional 
wind tunnel are discussed in reference 2. 

Since we will be interested in problems in which some of the nonlineari- 
ties may be large , conventional methods of analysis in which all nonlinearities 
are considered small (see, e.g., refs. 3 and 4) are not always suitable. Two 
methods that treat large nonlinearities in the static moment have been devised 
by Murphy: a perturbation technique (refs. 5 and 6) and a quasi-linear tech- 

nique (ref. 7) that yields results essentially the same as those of Coakley, 
Laitone, and Maas (ref. 8). Although these methods are quite general in that 
they are applicable for combined pitching and yawing motions , only cubic static 
moments are discussed, and explicit damping formulas for only linear darling in 
p lanar and near-circular oscillations are derived and plotted. Redd, Olsen, 
and Barton (ref. 9) derived formulas for the effective damping of nonlinear 
oscillations, but considered the static moment to be linear. 

Many aerodynamic systems are not describable by cubic static moments, and 
quintic or higher order polynomials are often essential for proper data reduc- 
tion. There is reason to expect that if the static moment requires description 
by higher order polynomials, then the damping moment will also. It is desir- 
able, therefore, to understand the influence of arbitrary nonlinear damping and 
static moments on oscillatory motions. The intent of the present paper is to 
contribute to this understanding. 

The analysis proceeds along two distinct lines after appropriate transfor- 
mations are used to render the original second-order differential equation 
into an equivalent integral equation on the one hand and an equivalent first- 
order differential equation on the other. By means of the integral equation, 
for which a special exact solution exists , an approximate solution is found 
that is valid for large damping nonlinearities and small effective static non- 
linearities. For small damping but arbitrarily large nonlinearities in the 
static moment, an approximate solution is then found by utilizing the first- 
order differential equation. For this latter line of analysis, the damping 
decrement can be conveniently represented by an "effective linear damping," 
which can be evaluated in terms of quadratures leading to elliptic integrals 
and gamma functions in special cases. Based on limiting values of the quadra- 
ture expressions, a simplified approximate formula for the effective linear 
damping is devised that is extremely accurate for a wide range of nonlinear 
parameters and, further, is amenable to data-reduction techniques. 


NOMENC IATURE 


A,B constants of integration 

a amplitude variable defined by equation (28b) 

b defined by equation (8lb) 

D 0 defined by equation (57c) 


2 


! 


,E congjlete elliptic integral of second kind 

f arbitrary function occurring in equation (27) 

H £ H s y s 

s>o 

Hpg effective linear damping defined by equation ( 66 ) 

H s dancing-moment coefficients 

Hg defined by equation (ll) 

h H 0 /v 

h s 2 + s y a (5*0) 

I defined by equation (8lc) 

K complete elliptic integral of first kind 

k modulus of elliptic integrals 

M ^M s y S+1 

s>° 

Me 

M s static-moment coefficients 

M s M s - H s 

m s M s y m S /v 2 or M s y a S /v 2 (eqs. (42), (54) ) 

P defined by equation (l 8 ) 

Q defined by equation ( 30 a) 

R defined by equation (34) 

s summation index 

t independent variable in equation (l) 

t Q initial value of t 

u,z functions defined by equations ( 3 ), ( 7 ), and ( 8 ) 


3 


y 

Y o 

Y m 

y 

y a 

y m 

y n 

yn+i 

y Q 

y (o) 

7 

9 

M- 

v 

T 

GJ 


y 2 /y a 2 

value of Y for zero damping 
rath order expansion term in Y (eq. (56)) 
dependent variable and amplitude of non! inear oscillation 
ii/2 


i (y/ + yJ+i)] 


a particular maximum amplitude 

nth maximum amplitude 

(n+l)th maximum amplitude 

initial amplitude of y 

zeroth order approximation (eq. (50)) 

defined by equation (79b) 

tan -1 (H 0 / 2 v) or sin' 1 (H 0 / 2 s /i^) 

m 2 (2m 0 + m 2 ) 

|Mo| l/2 if l^y 11 ! 172 if % = 0 

frequency variable defined by equation (28c) 
angular frequency 


BASIC ANALYSIS 


The Basic Equation 

In this study we wish to consider oscillations that are influenced by non- 
linearities in both the damping and static forces or moments. We shall assume 
that the nonlinearities are functions of the dependent variable only. 


Consider a nonlinear oscillation to be governed by the following differ- 
ential equation for y = y(t): 
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3?he middle term represents the damping moment, the last term, the static 
restoring moment. We shall refer to the term represented by s = 0 in both the 
damping and static moments as the linear moment. Similarly, we shall refer to 
polynomials terminated by s = 1, s = 2, and s = 4 as quadratic, cubic, and 
quintic moments, respectively. When H s and M s are nonzero only for even 
values of s, the motion is called symmetric; when H s and M s exist for odd 
values of s, the motion is called unsymmetric . Symmetric motions are of the 
greater interest, but unsymmetric motions, such as oscillations about a trim 
angle in aerodynamics, do occur. 

Equation (l) includes as special cases several well-known equations. If 
Hq^O, Mq^O, but H s = M s = 0 for s ^ 0, equation (l) is linear and has a 
well-known solution. If Hq = -H 2 = -e and H s = 0, s ^ 0, 2; Mq^O and 
M s =0, s / 0, then equation (l) reduces to the classical Van der Pol equation 

+ e(y 2 - 1) fjr + = 0 (2) 

If all the damping terms vanish , leaving only a nonlinear restoring moment, 
then the motion is periodic. If, in addition, the restoring moment, the last 
term in equation (l) , is a quintic or less in y, then exact solutions may he 
found in terms of elliptic functions for symmetric oscillations. 


Equivalent Pair of Equations 

It is convenient to seek a transformation that will express equation (l) 
in such a form that the first derivative of the dependent variable does not 
appear . 

Let us introduce two functions u(t) and z(t) so that 


y(t) = u(t)z(t) 


(3) 


Since we have introduced two functions to replace one, we have an arbitrary 
condition at our disposal. We will select u(t) so that the differential 
equation for z(t) does not explicitly involve its first derivative. 


If the transformation given by equation ( 3 ) is substituted into 
equation (l), then 


d 2 z 

dt 2 


s** y E 

udt L, 

v s>o 


U S S 
H S U Z 


dz 

dt 


a 2 u 

dt 2 


K 


ML 


du 

dt 


s>o 


¥ 


s s 
z 


z = 0 


( 4 ) 


It does not suffice to set the coefficient of dz/dt in equation (4) equal to 
zero in order to specify the needed equation for u. If this were done, a 
term involving dz/dt would appear in the coefficient of z owing to the 
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presence of d^/dt 2 . The form of the coefficient of dz/dt, however, suggests 
the following relation for u: 


1 du. 1 \ TT s s / _ \ 

u dt = " 2 /_ j a s H sU z (5) 

s>o 

Here the a s are constants to be determined so that the dz/dt terms vanish 
in equation (4) . 

Taking the derivative of equation ( 5 ) to determine d 2 u/dt 2 and substitut- 
ing into equation (4 ) , we find 





H s u s z s ; 


dz 

dt 



\ H s du 

__s s 

/u dt 

£ 

tSJ 




z = 0 


(6) 


We now eliminate the derivative dz/dt from equation (6) by requiring that its 
coefficient vanish. This can be done term by term if 


2 

a q — r 
s 2 + s 


Hence, 


d 2 z 

dt 2 



+ 


1 du\ 
2 + s u dt/ 


u s z s 


z 


= 0 


(7) 


1 du _ _ 
u dt 



s s 
u z 


s>o 


( 8 ) 


Equations (7) and (8) are the ones we have sought. These equations 
together with equation (3) replace the single equation (l) for y by two 
simultaneous equations in the new variables u and z. It may appear that the 
solution of the two simultaneous equations for u and z would be more compli- 
cated than the single equation for y. We will find, however, that the equa- 
tions become uncoupled for special values of M g and H s , and exact solutions 
may then be found. Moreover, transformations of these equations will yield a 
single integral equation for y and a single first-order equation for y that 
are appropriate for establishing approximate solutions. 


To see more clearly the nature of equations ( 7 ) and (8) , let us write 
them in a slightly different form. Substituting equation (8) into equation (7) 
gives 
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presence of d 2 u/dt 2 . The form of the coefficient of dz/dt, however, suggests 
the following relation for u: 

IdU 1 V x, S S ( c -\ 

u dt = ’ 2 /_ j a s H s u z (5) 

s>o 

Here the a s are constants to be determined so that the dz/dt terms vanish 
in equation (4) . 

Taking the derivative of equation ( 5 ) to determine d 2 u/dt 2 and substitut- 
ing into equation (4), we find 





H s u s z s 


dz 

dt 


+ 



\ H s du 

s si 

ju dt 

c 

tSJ 




z = 0 


( 6 ) 


We now eliminate the derivative dz/dt from equation (6) by requiring that its 
coefficient vanish. This can be done term by term if 


2 

a o — "r 
s 2 + s 


Hence, 


d 2 z 

dt 2 


+ 



+ 


1 du\ 
2 + s u dt / 


u s z s 


z 


= 0 


1 du _ 
u dt 


I 

s>o 


H 


2 + s 


s s 
u z 


(7) 

( 8 ) 


Equations (7) and ( 8 ) are the ones we have sought. These equations 
together with equation ( 3 ) replace the single equation (l) for y by two 
simultaneous equations in the new variables u and z. It may appear that the 
solution of the two simultaneous equations for u and z would be more compli- 
cated than the single equation for y. We will find, however, that the equa- 
tions become uncoupled for special values of M g and H s , and exact solutions 
may then be found. Moreover, transformations of these equations will yield a 
single integral equation for y and a single first-order equation for y that 
are appropriate for establishing approximate solutions. 


To see more clearly the nature of equations (7) and (8), let us write 
them in a slightly different form. Substituting equation (8) into equation ( 7 ) 
gives 
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i 




2 + s uSzS ) U - 0 


Equation (9) may be simplified by introducing new parameters, H„, defined by 


H s u s z s = 


H \ 2 

— u s z s \ 

2 + s 


The square of the series in equation (10) may be arranged so that H g may be 
identified as 


H s V Hnfe-n 
s L (2 + m)(2 + 


Equation (9) becomes 


^-§ + I ) (M r - Hju s z s z = 0 


With equation (T) written in this form, special exact solutions suggest 
themselves. 


Special Exact Solutions 


It can be seen that when 


M s = H s , s 4 0 


equation ( 12 ) reduces to 


<vV>- 


which has the general solution 


z = cos t 


where 


= v(t - to), V 2 = Mo 
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and t Q is an arbitrary constant. The other arbitrary constant which would » 
normally appear has been suppressed since it will appear in the solution of the 
equation for u(t), equation (8). 

Consider now a two-term damping moment given by 



ay _ 

dt 



(16) 


i 

Here we have kept the linear term H c and one nonlinear term of arbitrary | 

power n, where n is an integer. With this damping moment, equation (8) , 

takes the form 

1 du Hq 

u dt 2 

A solution to equation (l7) may be found by introducing the transformation 


Hn 


2 + n 


u n z n 


(17) 


= -H t 

r = P e ° 


(18) 


so that 


^IdP 2H n 

dt 2 + n 


(19) 


which has the solution 


H rt t 


p = U 2 e °~ = I A + 


r n 
2 + n J 0 


t -|H 0 t' 


2 

n 


dt 1 


( 20 ) 


where A is a constant of integration. 


From equations (10) and (13) we find that for this special case the only 
nonzero values of M s and H s are given by 

■'v 


= H„ = 


H 0 H n 


n “ n " 2 + n 


M an = H 2n = 


H, 


n 


(2 + n)‘ 


( 21 ) 


Finally, by combining the solutions for 
function y(t) can be expressed as a function 


u(t) and z(t), we find that the 
of t = v(t - t 0 ) by 
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y(t) = 


Hot 
2 V 


COS T 


ntt 


■n 


y Q 


(2 + n) v Jo 


r 

Jo 


Hq t 1 


cos 


dT' 


( 22 ) 


where y Q = y(t = t Q ) and t Q are arbitrary constants, and y(t) is a general 
solution to the following equation: 



(Ho + H ny n )|f + 


M_ 


< _HoHn n + H n 2 
^ 2 + n (2 + n ) 2 



(23) 


where n is an integer. 

Solution (22) to equation ( 23 ) "was first given by Smith (ref. 10) , 
although by an entirely different derivation. Subsequently, we shall use this 
exact solution as a basis of an approximate solution that does not restrict 
the nonlinear restoring terms Mq and M £n to the values given by equation (21). 

It is interesting to see the explicit forms that solution (22) takes for 
specific values of n. Here we present the solutions for n = 2 and n = 4, 
which are most likely to be of practical interest. 


n = 2 : 


y = 


. hi 

e 2V cos t 


H_ 


H q t 


B - e v [1 - sin 9 sin(2T - 9 )] 
kHo 


(2b) 


H ^ 

y2 = M o IT » T = v(t - t Q ) 

H 0 


sin 9 = 




n = 4: 


y = 


Hot 

2V 

e cos t 


B - 


H 


2Hot 


8Hq 


* e v 


1 - 


^ sin sin(2T - ty) 


i sin 9 sin(4r - 9 ) 


j. 

I 4 


( 25 ) 


sin ^ = 




■F + H c 


, sin cp = 


H. 
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where B is a constant. Formulas (24) and (25) are useful because they indi- 
cate explicitly and analytically how the motion may be affected by arbitrarily 
large damping nonlinearities. 

The special exact solution (22) has the shortcoming that the nonlinear 
static-moment terms ^ are not open constants , but are determined by the 
damping-moment terms H 0 and H n as given by equation (2l) . We would like to 
lift this restriction and find approximate solutions for arbitrary values of 
Mh- To do this, we must take into account that the frequency may be variable 
and recast the equations into more appropriate forms. 


Equivalent Integral Equation 


The 
which we 


two simultaneous equations that we shall deal with are (l2) and (8), 
rewrite here as 



(26a) 


1 du 2 
u 2 dt 


- 2 E 

s>o 


IL 


s^s 


2 + s 


trz 


(26b) 


where M s = M s - H s . The parameters M g may be viewed as the "effective" 
static-moment parameters. 


We introduce a modified version of the Kryloff-Bogoliuboff technique, 
which is outlined in appendix A. Consider a nonlinear differential equation 
having the form 

+ v 2 x + f(x, §) = 0 (27) 


where f is an arbitrary function of x and dx/dt, and v 2 is a positive 
quantity. If we assume that equation (27) has a solution of the form 


x = a(i-)cos t , t = T(t) 

then, as shown in appendix A, a(r) and x(t) are determined by 

f(a cos t, -av sin T)sin t 
da = _f_ 

dT , 1 j,/ . \ 

1 + — ^ f(a cos t, -av sin t;cos t 
av 


(28a) 


(28b) 
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(28c) 


* 


dr 

dt 


1 + — — f(a cos 


av £ 


t, 


-av sin t)cos t 


The amplitude variable, a, is determined as a function of t, the frequency- 
variable, by solving equation (28b). Using this solution, a = a(-r) , one may 
then find the relationship between t and the original independent variable 
t by solving equation (28c) in the form of a quadrature 


t 



dr* 



COS T* 


( 28 d) 


Consider now the application of the modified Kryloff-Bogoliuboff equations 
to our particular problem. To account for oscillations for which _Mq < 0, we 
shall select special values for the arbitrary parameter v 2 . If Mq £ 0, then 
we shal 1 choose v 2 = |Mq j . If Mq = 0, then we shall choose v 2 = |Mny m n |, 
where y m is a suitable characteristic amplitude of the motion. Noting that 
z = a cos r, we may identify the function f(a cos t, -av sin t) from equa- 
tion (26a) as follows: 


f = 




a s u s 



v 2 a cos t 


(29) 


Expression (29), in this case, is not in convenient form to solve equa- 
tion (28b) because the variable u is present, and u must be found as a 
function of a(r) and t by combining equations (26b) and (28c) to eliminate 
the variable t. Because this step cannot be easily accomplished, we instead 
introduce a new variable, Q, defined as follows: 

HqT 

Q = a^u 2 e V (30a) 


It now can be established that y is determined as a function of t by 

_ Hqt 

y 2 (r) = Q(r)e v cos 2 t (30b) 

An equation for Q(t) now remains to be found. Prom equation (30a) it 
follows that 


Hqt 


..2 —v~ da£ _ du 2 dt _ ^oQ 

dT dT u 2 dt dr v 


(30c) 


Substituting equations (26b), (28b), and (28c) into equation (30c), we obtain 
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_1_ dQ 
Q 2 dr 



This equation is the more general counterpart of equation (19) . The appropri- 
ate expression for f in terms of Q and t is found from equation (29) to he 



Thus, when equation (31h) is used for f in equation (31a), a first-order non- 
linear equation for Q(t) results. 

It is, in general, too difficult to integrate equation (31a) for 
Q = Q(t). Nevertheless, it is convenient to go a step further and find an 
integral equation for y that will serve as a basis for approximate solutions. 
We first note that f can be expressed in terms of y 2 by using 
equation (30b): 




where 


R(t) = — cos 4 t + 2 sin t cos 3 t (3*0 


Direct integration of equation (33) gives 



(35) 


Substitution of equation (35) into equation (30b) now yields an integral equa- 
tion for y(-r): 


Hot 



(36) 


Here we have identified the constant of integration B as l/y Q 2 , where 

y 0 = y(T = o) . 

If it is assumed that y = y(>r) can be found, then, at least formally, the 
corresponding relation for t(-r) can be determined from equation ( 28 d) : 



The integral equation ( 36 ) is interesting because it reflects roughly the 
nature of a solution that one might expect for a nonlinear oscillation- 
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Furthermore, when v 2 =M o >0, M s = 0 for s > 1, and H s = 0 for s ^ 0,2, the* 
exact solution (22) for a = 2 is recovered immediately. Equation (36) will 
thus be useful for a perturbation of this exact solution. 


Equivalent First-Order Equation 

We shall now derive the most useful and most important equation of this 
analysis. Since y = uz, the first derivative of y can be written in terms 
of u and z as 

d y - v/i . 1 /_os 

dt \z dt u dty ) 

We now again introduce the modified Kryloff-Bogoliuboff equations, derived in 
appendix A, so that z = a cos t and dz/dt = -av sin t- In addition, we make 
use of equations (26b), (28c), and (32b), and write equation (38) as follows: 



s>o 


The general second-order differential equation (l) with t as the inde- 
pendent variable has now been reduced to a single first-order equation with t 
as the independent variable . The relation between t and t is again provided 
by equation (37) • Although equation (39) is equivalent to the integral equa- 
tion (36) because it describes the same nonlinear oscillation, it has a very 
different appearance . 


Distortion Due to Nonlinear Damping 

Equation (39) can immediately be used to determine where the maximum 
amplitudes occur in a nonlinear oscillation. When there is no damping in an 
oscillation, the motion is periodic, and it is easy to deduce that the maximum 
amplitudes occur at T m = nit, where n is an integer. For a general oscilla- 
tion, the maximum amplitude occurs when dy/dT vanishes, and thus from equa- 
tion (39) we deduce that 

tan r m = - ^ 2^ (4o) 

s>o 


where y m = y(r m ) is the maximum amplitude. 



* Formula (40) indicates how much the maximum amplitude is shifted from the 
point of symmetry, T m = nit. Consider a special simple example, namely, the 
limit cycle for the exact solution (24). The limit cycle occurs when H Q < 0, 
H 2 > 0, and t -*• 00 and, in this case, equation (24) yields 


(y) 


lim 


/ 


COS T 


H. 


, — [1 - sin cp sin(2T - 9) J 
4Hq 


(4l) 


where 


or 


sin 9 = 


Hq 

2^ 


R 0 

tan q> = -2 
2v 


The value of y m from equation (4l) is 



and equation (4o) now yields for this example 


tan t, 


m 





= tan 9 


Hence , 

t = rut + cp 
m 



Figure 1.- Distortion of motion due to 
nonlinear damping. 


and the limit motion is distorted from 
the symmetry position by an amount 
9 = tan _1 (H 0 /2v) . This motion is shown 
in figure 1 for values of 9=0°, 

-15°, -45°, and -75°. When |9I is 
large, a substantial distortion of the 
nearly cosine oscillation characteris- 
tic of small nonlinearities is found, 
and the oscillation is no longer sym- 
metric about T m = nit. 
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EXACT SOLUTIONS FOR ZERO DAMPING 


It is interesting and instructive to see how exact solutions for zero 
damping appear in the present formulation of nonlinear oscillations. Moreover, 
these exact solutions will be useful later in approximations for nonzero 
damp ing . 


General Solution 

Let us now consider only symmetric motions, so that s is even, and all 
the damping terms are zero (H s = 0 ) . The motion will be periodic and have a 
maximum amplitude that we denote by y m - We introduce a normalized form of y 
by the transformation 



In addition, we introduce the coefficients m s defined by: 


_ _ M sy m s _ > 0 

m s - ' s ^ 0 


( 42 ) 


Note that if M Q / 0, m g is a ratio of the magnitude of the sth nonlinear 
term of the static moment to the magnitude of the linear term. With this nota- 
tion, one may express equation (39) for zero damping as follows: 



-2Y 0 tan t 


. 2 

sin 


T + COS 2 



m s^o 


S/2 


s>o 


(43) 


Equation (43) can be immediately integrated if a transformation to a new 
independent variable x = tan t is made. If the boundary condition is imposed 
that Y 0 = 1 when t = nrt (n is an integer), then 


tan 2 



(44) 


The zero-damping solution becomes complete when a relationship between t 
and t is established. This is done by means of equation (37) j which now 
becomes 


l6 



T* + COS 2 T* 


v S/2 
m s Y o 


Formula (45) is not always a convenient expression to use because Y Q is gen- 
erally an implicit function of r. Consequently, an alternate expression in 
which t is determined as a function of Y Q is often more useful. The alter- 
nate formula may be found by changing variables according to equations (43) and 
(44) , and thus obtaining 


t - t„^ = 


2 + s 


(46a) 


nfl < T < nfl + g and 0 < Y q < 1 

Expression (46a) has been obtained by more conventional methods. If 
m s = 0 when s > 6, then equation (46a) can be evaluated by elliptic integrals 
(see refs. 11 and 12). The angular frequency w of the oscillation is defined 

■by 


v _ 1 

o> Jt 


2 + s 


(46b) 


An approximate formula for the frequency with arbitrary values of m s has been 
derived in reference 13 * 


Special Case of a Cubic Static Moment 

For the case of a linear -plus- cubic moment, m s = 0 when s > 4, it is 
easy to find an explicit function for Y 0 (t) from equation (44); namely. 


IT 


(p/m 2 )cos 2 * 


T 


(47)* 


Y o " 


sin 


t + m 0 cos 


T + J{ 


sin 2 t + m 0 cos 2 t) 


+ p cos T 


where p = m 2 (2mQ + ni 2 ) • 

With this explicit representation, equation (37) yields the following 
relation between t and t 


t - *0 " T 


d T ' 


'o 


Ji sin 2 t' + m 0 cos 2 t ' ) 2 + p cos 


(48a) 


4 t’ 


It follows that the frequency is given by 


v 

0) 



dr 


J( sin 2 t + m 0 cos 2 t ) 2 + p cos 4 t 


(48b) 


It is interesting to compare equation (48b) for the frequency with the 
corresponding formula (46b), which may be evaluated in terms of elliptic inte- 
grals. The quadrature (48b) is a proper integral and may be evaluated 
numerically without difficulty. 


Formulas (48) contain four possible cases of interest, depending on the 
value of the moment terms Mq and M 2 . To discuss them, we adopt the following 
convention. If Mjj is positive, it is called stable; if it is negative, it 
is called unstable. _The various possibilities are stated by describing the 
term Hq and then M 2 . For a cubic static moment we thus have the following 
cases : 2 


Stable-stable and stable -unstable static moment 


M q > 0, v 2 = M q , m 0 = 1 

Unstable -stable static moment 
Mq < 0, v 2 = |Mq| , m Q = -1 
Pure cubic static moment 

Mq = 0, v2 = ^m 2 > 0, m 0 = 0, m 2 = 1 


0 < m 2 < oo (stable-stable) 
-1 < m 2 < 0 (stable -‘unstable) 




2 < m 2 < 


> (49) 


2 Note that if m 2 < -1 for the stable-unstable case, tumbling motion 

results; also, if 1 < m 2 < 2 for the unstable -stable case, unsymmetric 
motion about a trim angle occurs. Since neither of these cases is covered in 

the present analysis, initial conditions must be such that the restrictions on 
m 2 are .observed. 
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The oscillations (eq. (47)) for the ciibic moments are plotted in figures 2(a), 
*(b), a nrl (c) to show the effect of m 2 for the various cases. In these fig- 
ures, n /Iq = y/y m is plotted versus the transformed independent variable, t, 
over a q uar ter cycle of oscillation. As indicated by equation (42), m 2 is a 
ratio of the magnitude of the nonlinear part of the static moment to the magni- 
tude of the linear part. The value m 2 = 0 represents the pure linear case, 
m 2 = oo the pure cubic case. In figure 2(d), y/y m is plotted versus the orig- 
inal independent variable, t, normalized in such a way that the curves have the 
game, frequency. This is equivalent to comparing each nonlinear moment with a 
different linear moment, in each case the linear moment being chosen to give 
the same frequency of oscillation as did the particular nonlinear moment . With 





(c) Variation with t; 

unstable -stable 1-3 

(d) Variation with t 

static 

moment . 



Figure 2*- Effect of 

on oscillatory 


behavior . 




this kind of plot, it can be seen that most nonlinearities in the moment have 
no drastic effect on the appearance of the motion (although the frequency is 
greatly affected). For instance, all stable-stable cases would fall within the 
shaded region of figure 2(d). Furthermore, the vast majority of stable- 
unstable cases ever encountered would be bracketed by the curves labeled 0.0 
and -0.8; most unstable- stable cases would be bracketed by the curves labeled 
2.5 and 00 . 


APPROXIMATE SOLUTIONS 


We now wish to find solutions for arbitrary values of H g and M s . With 
the use of the special exact solutions already discussed, it is possible to 
find approximate solutions valid when the parameters H s and M s are small in 
a certain sense. Essentially, we shall find first-order perturbation solu- 
tions to the exact solutions that are special cases of equations (36) and (39) • 
We shall consider only symmetric motion and hence make the restriction 
H s = M s = 0 when s is odd. 


Perturbation Solution Using the Integral Equation 

Because the integral equation (36) has an exact solution for a special 
case, that is, when Mq > 0 and M s =0, s>0, and when H s = 0, s > 4, it is 
natural to ask if there exists a small perturbation of this exact solution in 
which certain parameters are small in an appropriate sense. The exact solu- 
tion can be expressed in the form of equation (22), with n = 2, 



Hq t 


v 

e cos 



T 


cos 2 


d T ' 


(50) 


Expression (50) can be considered as a zeroth order approximation. Since equa- 
tion (36) is an integral equation, it is straightforward to use equation (50) 
and determine a first-order approximation to equation (36) as a first step in a 
scheme of successive approximations. 

Consistent with the zeroth approximation, we require that M Q > 0, and 

that 


and 


1 

S>2 



« 1 


(51a) 



^H s 

(2 + s)H 2 


S-2 


y 


« 1 


(51b) 
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.With these conditions, we can expand the denominator of the integrand of the 
integral in equation (36) in a Taylor series and keep the first term- The 
first-order approximation can now he obtained by replacing the y 2 that 
appears under the integral by the zeroth-order approximation y( Q ) and retain- 
ing terms consistent with equations (5l)- The first-order approximation thus 
can be written as 



(52) 


Formula (52) for the first approximation is expressed in terms of quadra- 
tures involving the zeroth-order approximation. In general, these quadratures 
are difficult to evaluate in terms of tabulated functions. For the special 
case when the nonlinear damping terms are zero, however, the quadratures can be 
evaluated with the use of elementary functions. We now write the corresponding 
formal solution for the frequency function t as a function of t. Starting 
with equation ( 37 ) and using the approximations consistent with the first- 
order solution, we obtain 


v(t - t 0 ) 



dr 1 


(53a) 


where, as before, y| Q j is given by equation (50) . Since higher order products 

of M g are neglected in the first approximation, we can also express t 
explicitly in terms of t by rearranging equation (53a) to read 


T = 



(53b) 


In many practical problems of interest in aerodynamics, the nonlinear 
re storing -moment terms are not small, and thus the above first-order approxima 
tions are not particularly useful. The damping, however, is usually small; 
consequently, results having wide application can be obtained from a perturba- 
tion solution of the first-order differential equation ( 39 ) • 
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Perturbation Solution Using the First-Order 
Differential Equation 


Let us write equation (39) in a normalized form that naturally suggests a 
perturbation procedure . We introduce a characteristic amplitude y a that we 
shall identify later and define the following symbols: 


Y = 



h = 



v 


m 


s 



h 


s 


1 H s S 

— ii; y a 


s even 


( 5*0 


Since M g « M g , m s is essentially a ratio of the magnitude of the sth non- 
linear term of the static moment to the magnitude of the linear term; h g is a 
similar ratio for the damping moment but multiplied by a constant, l/(2 + s) . 
With this notation, equation (39) becomes 



(55) 


Since an exact solution exists for equation (55) when h = 0, namely, 
equation (44), we seek a perturbation solution of equation (55) valid when 
h « 1. Therefore, we consider that h « 1 and h g = 0(l) and assume over a 
half cycle, at least, that Y(t) can be expanded in a series of the following 
form: 


00 



m=o 


(56) 


Substituting equation (56) into (55) > collecting terms, and requiring that the 
coefficient of each power of h vanish identically yields an infinite set of 
equations for the functions Y m . We shall confine our interest to the first- 
order perturbation and thus to only the functions Y Q and Y x . The equations 
for Y 0 and Y x are found to be 


dY 0 _ 2Y 0 tan t 

d T D q 


(57a) 
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. 2 tan t 
D ° dr + 


o 


sin 2 


t + cos^ T 


Z(-t) -.c i x - -a,, I vr 

s>o J s>o 


where 


D q = sin 2 r + cos 2 t 


I “sC 


(57b) 


(57c) 


s>o 


With the use of equations (57a) and (57c), it is not difficult to verify that 
an integrating factor for equation ( 57 b) is merely sec 2 t, so that equa- 
tion ( 57 b) can he expressed as 


■A(D 0 sec 2 r Y x ) 


-2 


sec 2 t Y 0 


I 

s>o 


* v s/2 
h s Y o 


(58) 


With the above equations we will be interested in the solution between two 
consecutive maximum amplitudes, say the nth maximum amplitude, y , and the 
(n+i)th, y n+1 * The zeroth order equation ( 57 a) yields a periodic solution 
with a given maximum amplitude. We now select y a so that the maximum ampli- 
tude squared of the zeroth order solution is the mean of the squares of the 
nth and (n+i)th maximum amplitudes 


y a 2 = ipi 



+ yit+i) 


We can now take equation (44) to be the 
solution of ( 57 a) and identify y m 
with y a • The situation is shown in 
sketch (a) . 

Correct within the first-order 
approximation, y n occurs at t = nfl 
and y n+1 occurs at t = (n + l)«. 

We integrate equation ( 58 ) so that 
Yjtnfl) = -Y^nfl + jt) 


Yx(r) = - 


2 COS 2 T 


r 

-Vi)* 


sec 2 t' 


••I 


h s Yo /2 dr’ 


s>o 


(59) 
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The function Y^t) is expressed by a quadrature. Since Y q (t) can be written- 
explicitly as a function of t only for special cases, such as the cubic 
restoring moment (eq. (47)), a general explicit form for the integrand cannot 
be found. With a change in variable, however, the quadrature can be expressed 
in the form of elliptic or hyperelliptic integrals. 


Consider a change of variable in the integral of equation (59) to Y 0 , as 
defined by equations (57a) and (44). Then, in terms of Y Q , the function Y^^ 
can be written as 



for nit < t < (n + l)it. The plus sign is used when [n + (l/2)]it < t < (n + l)it 
and the minus sign when nit < t < [n + (l/2)]it. It is clear that Y x (t) and 

G(Y 0 ) are odd functions of {t - [n + (l/2)]it}. When the polynomial in the 
radical of equation (60b) is a quartic or less, then the integral can be eval- 
uated in terms of elliptic integrals and elliptic functions. We shall eval- 
uate explicitly the complete integral later for special cases of the damping 
decrement . 

The corresponding first-order form of the frequency equation (37) can now 
be determined. Inserting equation (56) into (37) and keeping only first-order 
terms in h yields 


v(t - 



dT* 


sin 2 t* + cos 2 t 1 ) fli s Y 0 


S/2 


s>o 


S-2 


COS 2 T 1 Y-J_ . — s *o 

s>6 


sm^Y^ dT 1 


+ 0(h 2 ) 


sin 2 t 1 + cos 2 t 1 


(n+^)it 


s>o 


(61) 
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where here t = t Q when t = [n + (l/2)]fl. This completes the formal first- 
order solution. Because the first-order integral in equation (6l) is very com- 
plicated , we shall not investigate it further here , but merely note that it can 
be integrated numerically without great difficulty when Y 0 and Y x are known 
functions of t . 


THE DAMPING DECREMENT 


We now wish to derive formulas that portray the distinguishing features of 
the nonlinear oscillation. Of an oscillation that is sufficiently well-defined, 
two distinguishing features are the maximum amplitudes and the half period. 

(By half period, we mean the time interval between two successive maximum ampli- 
tudes of the motion. In this sense, it is a pseudo-period, since the motion is 
not necessarily periodic. Given the half period, we can likewise define an 
angular frequency.) A pertinent parameter to consider is the ratio of two suc- 
cessive maximum amplitudes, called the damping decrement. A slightly more sig- 
nificant parameter is the logarithmic decrement because it is a constant for a 
strictly linear oscillation and thus forms a good standard for comparison. 


Decrement for Small Damping 

Consider now the first-order solution that neglects terms of order h 2 in 
equation (56); that is, 

Y(t) = Y 0 (t) + hY-,0) + 0(h 2 ) (62a) 


From this expression we find y n 



and y n+1 to be given by 
1 + hY 1 (mt) + 0(h 2 ) 

1 + hY^nit + it) + 0(h 2 ) 


(62b) 

(62c) 


Noting that Y x (nit) = -Y 1 (nit + it), dividing one equation by the other, and 
expanding for small h, we find the damping decrement to be 



1 - 2hY 1 (nit + it) + 0(h 2 ) 


(63) 


Taking the logarithm of this expression and again expanding for small h 
yields the following formula for the logarithmic decrement: 


Zn 


y n 

y n +i 


= -hY 1 (nit + it) + 0(h 2 ) 


( 64 ) 


25 


Here Y^ntf + Jt) can be evaluated by equation (59) or (60) . Before proceeding,, 
however, it is convenient to modify equation (64) by defining an effective 
linear damping. 


Effective Linear Damping 


In a strictly linear oscillation it is known that the logarithmic decre- 
ment is a constant with the value 


l n 


y n 

Yn+i 


Ho? 

2v 


(65) 


where v is the frequency defined by v 2 
an effective linear damping by 


Zn 


yn 

y n +i 


_ H °e* 
“ 2 w 



Analogously, we can define 


V it v 

h 2 « 



( 66 ) 


where now, consistent with the first-order analysis, we define the frequency 
w by 


or 



(67a) 


(67b) 


^The logarithmic decrement can be written in more conventional form if we 
introduce the so-called damping ratio & as & = H 0 /2M 0 . Then, realizing that 

v = = A, - <h 0 2 a> for a linear oscillator, we have 

t _ y n _ 7t6 

In - = ■ ■ 

yn+1 Ji - 6 2 

For a linear system, the logarithmic decrement is independent of amplitude, 
and, therefore, we can also write 


Zni 


^n 

Yn+2 


2 flS 

Jl - 6 2 
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- We now combine equations (6*0 and (66) to get 

= - — — Y (n« + it) + 0(h) (68) 

H It V 1 

o 

Expression (68) is an important formula, and we shall investigate it in some 
detail . 


General Formula for Effective Linear Damping 


The formal first-order value for the effective linear damping can be 
found by evaluating Y (n« + n) in formula (68) by (60) (noting that 


Y 0 (nit + «) = l) 


Ho 


2 w 
jt v 


I 

s>o 



Jo 


2 + sY 


T)l v; 

/ s>o 


S/2 


dY, 




(69) 


where to/v is given by equation (67b) • Formula ( 69 ) can be explicitly evalu- 
ated in terms of tabulated functions for restricted combinations of m s . We 
shall investigate these special solutions in order to gain an appreciation for 
the influence of the nonlinearities on the damping, and to ascertain the asyrmp 
totic values of the general formula in various limits. On the basis of this 
information, an approximate , but quite general, formula can be established 
that is valid for almost the whole range of possible values of the nonlinear 
parameters. 


Cubic Damping and Cubic Static Moment 


Formula ( 69 ) can be evaluated in terms of elliptic integrals for a cubic 
static moment; that is, m s = 0 when s > 4. For the sake of simplicity, we 
also consider only a cubic damping moment; that is, H g = 0 when s > 4. For 
this case, formula ( 69 ) appears as 




1<») f |m 0 + -£■ ^1 + Y Q ^ (l + 2h 2 Y D ) dY 0 


H q mo + m 2 


(70a) 


f Y 0 (l - Y 0 ) 


m 2 / \ 
mo + -jrd + Y q ) 


where 
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V 

oo 


1 

Jt 



1 



dY 


o 




(1 + Y q ) 


(70b) 


Two representations are required for the various static moments defined 
by equation ( 4 - 9 ) to express ( 70 a) in terms of elliptic integrals. (No such 
distinction is needed for the various cubic damping moments that could exist.) 
These are as follows: 


Ho e 2m 0 + m 2 


H q m Q + m 2 


A (x + 1 

2 V 2mo + m 2 K 


f +»B f + 


m 2 7. 


K 2m n + m 2 K 


(71a) 


K = K(k) = 


Jt/2 


dcp 


Jl - k 2 sin 2 cp 


Complete elliptic integral 
of first kind 


(71b) 


E = E(k) 




k 2 sin 2 9 <19 , 


Complete elliptic integral 
of second kind 


( 7 ic) 


Stable-stable and unstable- stable static moment: (formula 236. l6, ref. 12 ) 


7, = 75 t(k 2 - 1)K + E] 

1 k 

7 = [( 3 k 4 - 5 k 2 + 2 )K + 2 ( 2 k 2 - l)E] 

2 3 k 

7 = -I5 [4(2k 2 - 1)7 + 3(1 - k 2 ) 7 ] 

3 5 k 2 1 

2 _ 1 

K 2 m 0 + m 2 



Stable-unstable static moment: (formula 23^.16, ref. 12) 


28 


7 X = h (K - E) 

k 

7 = “^4 [(2+ k 2 )K - 2(1 + k 2 )E] 

2 3k 

7 3 = t^(k 2 + 1)7 2 - 37 J 

5k 

k = 

2 + m 2 

An equivalent expression for H 0 e /H 0 for a cubic static moment, but a 
linear dancing moment, that is, h 2 = 0, has been given by Murphy (ref. 5) and 
Murphy and Hodes (ref. 6). Murphy obtained his results in a different manner . 
He assumed that the damping was small and then used the exact elliptic func- 
tions obtained for zero damping in an averaging process over a cycle of motion. 
The perturbation method of Murphy thus yields the same value for the effective 
danping as the present analysis for a cubic static moment, but it is difficult 
to extend the method of Murphy to higher order static moments owing to the 
complexities in representing more complicated static-moment solutions by ellip- 
tic functions. An approximate formula for Ho e /H 0 for a cubic static moment 
and linear danping has also been derived by Murphy (ref. f) by a ”quasi-linear 11 
method, and we will discuss this later in the section entitled, "Comparison for 
a Cubic Moment.” 



Expression (7la) for H 0e /H 0 in terms of elliptic integrals is actually 
a fairly involved formula. Also, the quadratures given by equation (70) are 
difficult to evaluate numerically because of their singular nature at both 
limits. It is interesting that the effective linear damping can be expressed 
by another set of quadratures that possibly may be easier to evaluate numeri- 
cally (especially with electronic computers) than either equation (70a) or 
(7la) . We obtain this alternate formula by transforming equation (70a) from 
the Y q variables to the t variables with the use of equations (43), (44), 
and (47) - The alternate result is 



(7^a) 


where 
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dT 


(74b) • 


V 

U 


i r 

°o 


Ji 


sin^ t + m n cos' 1 



+ l-i cos 4 T 


and 


li = m 2 (2m 0 + m 2 ) 


The above quadratures are proper integrals and can be evaluated numerically 
without difficulty. It is worth mentioning that equation (74a) can easily be 
extended to higher order damping terms. 

There are several ways in which formulas (71a) and (74a) can be plotted. 

If one is interested in knowing H 0e /H 0 for given values of h 2 and m 2 , carpet 
plots are useful. Such carpet plots are shown in figures 3(a), (b) , and (c). 
Because of the ease of making linear interpolations, these plots are convenient 
to use . 



(c) Unstable -stable 1-3 static moment* 


Figure 3*- Carpet plot of effective linear 
damping. 


30 


Another method of plotting illustrates how Hq /H q varies with amplitude 
for given values of the damping- and static-moment parameters. We thus plot 
HQ e /H 0 as a function of m 2 and set h 2 = (h 2 /m 2 )m 2 so that 

h 2 _ H 2 y 2 _ H 2 |Mq| 
m2 M 2 4h q m 2 

can he varied as a parameter that does not depend on amplitude for a given 
oscillator. These plots are illustrated in figures 4(a), (h), and (c). The 
curves denoted by h 2 /m 2 = 0 were given by Murphy and Hodes (ref. 6 ). The 







curve^s in figures 4(b) and (c) are singular for the stable -unstable moment at 
m 2 s M 2 y n 2 /Ko = "1 and for the unstable- stable moment at m 2 = M 2 y n 2 /|M 0 | = 2. 
The dashed lines are the asymptotes for m 2 «>• It is possible to show by an 
expansion for large m 2 that the asymptotes are given by 


where 


or 





m o^ 1 2 

m 2 



m 2 



f = gamma function 



0.667 + 


0 


296 


^oV\ 
m 2 J 


+ 0.73U 



(75) 


Likewise, for the stable-stable and stable-unstable moments, the variation as 
m 2 0 is given by 

“5” = 1 + (" 16 + mi) m2 + 0 ( ms2 } ( 76 ) 

Another interesting aspect of these curves is the appearance of maxima and 
minima in some cases. If these curves were established by experiment, this 
behavior might be considered perplexing without theoretical background. They 
portray the effects of strong interactions of the nonlinear damp ing and static 
moments. Still another interesting feature is the curve labeled h 2 /m 2 = 5/4 
in figure 4(b) . Detailed examination of equation (74) shows that this value 
is the dividing point between cases where a limit cycle can or cannot exist, 
given a stable-unstable static moment. For h 2 /m 2 < 5/4, there is no 
possibility of a limit cycle. 

All curves shown in figure 4 were compared with results obtained from 
numerical integrations of equation (l)^ In these numerical integrations, var- 
ious values of the parameter h = H 0 /|M 0 | 1/2 were chosen up to a value of 
about 0.1. This encompasses the practical range of damping encountered in 
ballistic-range testing. In all cases, the results of the numerical integra- 
tions agreed almost exactly with the results from formula (7la) or (74a) shown 
in figure 4. To the scale these curves are plotted, no differences can be 
shown . 

When H 0 < 0, a limit cycle exists when H Qe = 0. The amplitude of the 
limit cycle may be found by setting H Qe = 0 in equation (7la) or (74a). 

Using equation (74a) we get 
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(sin 2 t-hhq cos 2 r )+ J ( sio ^ T+m 0 cos 2 t) 2 +p cos 4 t 


2(2m 0 -ttn 2 ) f */ 2 
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Jo 


COS 2 T dT 


£(sin 2 t-wiq cos 2 


*r)-^/(sin 2 t-hoq cos 2 t) 2 +H cos 4 tJ 



( 77 ) 

An alternate formula is obtained from 
equation (7la), which is equivalent to 
the result found by Murphy and Hodes 
(ref • 6 ) by investigating the singu- 
lar nature of the limit cycle in an 
simp lit ude plane . The effect of the 
cubic static moment on the amplitude 
of the limit cycle is demonstrated in 
figure 5 * Except in the region of the 
singularities, the effect of m 2 is 
not large, yielding a value of 
(-h 2 )2 = 0.9118 as m 2 -*• °° contrasted 
to the well-known value (-h 2 ) 1 = 1 
when m 2 -*• 0 . It will be shown later 
that for quintic and higher order 
static moments, the effects are much 
greater . 


Arbitrary One-Term Static Moment 


The quadratures in equation (69) can be evaluated in terms of tabulated 
functions for the special case of an arbitrary one-term static moment . This 
case yields the asymptotic values that are obtained as the highest order term 
in an arbitrary polynomial static moment tends to large values. 


We now consider 
we choose v 2 = M n y n 
now appears as 



Ho 


1 

x 


M ^ 0 for s = n, where n is an even integer. Hence 
> S 0 so that m s = 0 , s ^ n, and m n = 1 . Expression (69) 
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and 


1/2 
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(78b) 


v _ 1 / 2 + n \ 

w m v 2 J 


dY 



n+ 2 ' 



o 


( n.+2 ) / 2 

With the substitution x = Yq , these quadratures can be written in the 

form of beta functions, which in turn can be evaluated in terms of gamma 
functions, T. For equation (78a) , we thus obtain 


where 



V 7 b 

/ , / sn s 
s >o 



S>2 


7 sn^s 


(79a) 


7sn 


4(2 + s) 
2s + 4 + n 


r(LtiWli_ + I s ) 
V * 2j \ H + 2 27 


s, n > 0 


(79b) 


Although formula (79) was derived for a one-term static moment, it also 
represents the asymptotic value of H 0 ^/H 0 for a polynomial static moment in 
which the nth term, m n > is dominant and stable and tends to infinity. For 
purely linear damping it is interesting that the asymptotic value is simply 
H 0e /H Q = 4/(n + 4). On the other hand, for a linear static moment, n = 0, 
equation (79a) can be written 


^ = i + 2 y — 

H o 1 2*4*6 


• (s - 1) 


S>2 


h 


(79c) 


Redd et al. (ref. ^) , considering only a linear static moment, computed 
the same value as equation (79c) up to s = 8 by a different method. Their 
conclusion, however, that a nonlinear static moment does not significantly 
alter the formula (79c), except near unstable trim, is not substantiated by the 
present analysis. Here it is found that equation (79c) is a special case of a 
more general formula (79a ) y and each type of nonlinear static moment may 
indeed affect H 0e /H 0 significantly. 

Expression (79a) may be used to examine the effects of the nonlinear 
static moment on the value of the limit amplitude. Consider a two-term damp- 
ing moment given by H s = 0 when s ^ 0 or k (k even) so that the first term 
is linear and the second term is nonlinear of arbitrary order. The limit cycle 
occurs when Hq < 0 and may be determined by setting H 0e = 0 to obtain 
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The value of 
values of (-h k ) j 


is a function of both k and n, and some of the 
are given in the following table: 


X k 2 

4 

6 

8 

0 

1.000 

1*333 

1.600 

1.829 

2 

.912 

1.167 

1.368 

1.540 

k 

.870 

1.084 

I.250 

1*392 

6 

.845 

1.034 

1.178 

1.300 

8 

.829 

1.001 

1.130 

1.238 

00 

.750 

.833 

.875 

.900 


It can be seen from this table that the value of (-1%) 1 decreases for a 
given value of k as n increases. Here we recall that n = 0 represents a 
purely linear static moment . 


Expansion for Small. Static-Moment Nonlinearities 

Consider the value of Ho e /Ho when the nonlinearities in the static 
moment are small. Here we assume M 0 > 0 so that iUq = 1, and we assume that 
m s is small so that higher order products and cross-products can be neglected. 
With these conditions, it is possible to expand equation (69) for small m s 
and evaluate each integral term by term. The result can be written as follows 


^e 

Ho 


I 

k>o 


Hk 


2l(k - l) 


-x 

S>2 


8ks m s 


+ higher order products of m c 


(8la) 


where 

b Vo = — - — [ (k + 2)l(s + k + l) + (s - k)l(k - l) - (s + 2)l(k - l)l(s + l)] 
ks 2 + s 

(8lb) 


l(n) 



,n/2 


dz 


(1 - z) 


1/2 


= 1 , 


n = -1 


1 * 3-9 » « « n 
2 * 4*6 . . . (n + l) 


n odd 


(8lc) 
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Some of the values of b ks for various values of k and s are tabulated as 
follows : 


X k 
sX 

0 

2 

4 

6 

8 

2 


1 

. , 21 

l 

109 


8 

4 

128 

8 

1024 

4 

11 

21 

_22 

12 

129 


12 

96 

128 

96 

768 

6 

279 

121 

697 

127 

3-3419 


256 

256 

2048 

512 

65536 

8 

193 

171 

_222 

1169 

77^2 


160 

320 

2560 

4096 

32768 

CO 

2 

1 

2 

2 

21 




4 

8 
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Expression ( 76 ) for the cubic damping and cubic static moment is a special 
case of ( 8 la) . Formula ( 8 la) will be of use in the construction of a general 
approximate formula for H 0 e /H Q . 


APPROXIMATE GENERAL FORMULA FOR H^/Hq 


Although the formal representation (69) is a very general result, it is 
difficult to evaluate for arbitrary combinations of the damping and static 
moments • Even for the cubic case, the evaluation in terms of elliptic inte- 
grals is rather complicated* Based on the special cases ( 79 ) and (8l), how- 
ever, it is possible to construct an approximate formula that is extremely 
accurate . 


A perusal of the formal solution (69) indicates that H 0e /H 0 
expressed by 



s>o 


can be 


(82) 


where a^ s and are functions of the parameters m 0 ,m 2 ^m 4 ,. . . . The 

functions a^ s and P^ s depend only weakly on these parameters, however, and 
thus we wish to find approximate values of c% s and P^s that are independent 
of m 0 ,m 2 ,m 4 ,. ... We can do this by requiring that equation (82) have the 
proper behavior when the nonlinear static moment is very large and stable and 
when the nonlinear static moment is very small; that is, equation (82) should 
have the same behavior as (79) and (8l) in these particular limits. 
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Formula (82) is "thus guaranteed to be valid when the dominant term 
m 0 ,m 2 ,m 4J , . . . 00 and also when the terms mo,m 2 ,m 4 , . . . -* 0. 


We begin, somewhat arbitrarily, by requiring that P^q = 1 - We do this 
in order to agree with the leading term of the sum ^_m s that appears in the 

denominator of equation (69). It can be shown that if equation (82) is to 
reduce to both (79) and (8l) in these respective limits, then a^ s and P^g 
are given by 


“to = 21 < k - 

_ ^ks^ks 
' “ko - 7 ks 


PkO - 1 




^ks 


ks 


“ko ” ^ks 


s > 2 [ 


( 83 ) 


where 7]j S > and l(n) are defined by equations ( T9t> ) , (8lb), and (8lc). 

Some values of a^ g and P^-g are presented in the following tables. 4 


“ks 



0 

1 

2 

3 

4 

5 

6 

7 

8 

0 

2.0000 

1.2732 

1.0000 

0.8488 

0.7500 

0.6791 

0.6250 

0.5821 

0.5469 

1 

1-53^9 

1.0095 

.8065 

.6913 

.6143 

.5582 

.51*17 

.4797 

.4510 

2 

1 .2500 

.8398 

.6797 

.5875 

.5250 

.4788 

.4428 

.4136 

.3895 

3 

I.O568 

.7208 

.5894 

.5130 

.4608 

.4218 

•3912 

.3663 

•3*+53 

4 

.9167 

.6323 

.5213 

.4565 

•4n8 

.3784 

.3519 

• 3301 

.3118 

5 

.8102 

*5639 

.4679 

.4n8 

.3731 

.3439 

.3206 

•3014 

.2851 

6 

.1266 

.5092 

.4249 

.3756 

.3*0-5 

.3156 

.2950 

.2779 

.2633 

7 

.6 589 

.4645 

.389^ 

.3456 

.3152 

.2920 

.2735 

.2581 

.2450 

8 

.6031 

.4272 

.3596 

.3202 

.2928 

.2720 

•2552 

.2413 

.2294 

OO 

0 

0 

0 

0 

0 

0 

0 

0 

0 


Note that a ks and P^ s are given for both even and odd values of k 
and s. We are still considering only symmetric oscillations, but can extend 
the analysis to cases where k and/or s are odd by insisting that both the 
static moment and damping moment be odd functions of the dependent variable • 
This is done by employing absolute values where needed. For instance, a qua- 
dratic static moment would be written as M = IV^y + M^y |y | • 
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^ks 


\k 
s 1 

0 

1 

2 

3 

4 

5 

6 

7 

8 

0 

1.0000 

1.0000 

1.0000 

1.0000 

1.0000 

1.0000 

1.0000 

1.0000 

1.0000 

1 

•9593 

.9575 

•95^2 

.9505 

.9468 

♦9 k 3 k 

•9399 

.9364 

.9337 

2 

•9375 

•93^6 

.9297 

.9243 

.9188 

•9135 

.9084 

.9038 

.8997 

3 

.9247 

.9211 

.9155 

.9091 

.9027 

.8964 

.8905 

.8850 

.8797 

4 

.9167 

.9127 

.9067 

.8998 

.8929 

.8861 

.8797 

•8735 

.8677 

5 

•9115 

.9073 

.9011 

.8940 

.8868 

.8798 

.8729 

.866 5 

.8603 

6 

.9082 

•9039 

.8975 

.8904 

.8831 

.8758 

.8688 

.8621 

.8558 

7 

.9060 

.9016 

.8953 

.8882 

.8808 

.8735 

.8664 

.8596 

.8530 

8 

.9047 

.9003 

.8940 

.8868 

.8795 

.8722 

.8650 

.8582 

.8516 

00 

1.0 

1.0 

1.0 

1.0 

1.0 

1.0 

1.0 

1.0 

1.0 


Comparison for a Cubic Moment 
For a cubic moment, equation (82) looks as follows 

Ho e _ mo + 0. 6250m 2 m Q + 0.6797m 2 /QI . 

— + ho + • • « ( 04 ) 

Ho m Q + 0.9375m 2 m Q + 0.9297m 2 

The accuracy of this expression can be demonstrated in several ways. First, 
we will compare equation (84) with the exact solution (7l) or (74) when only 
linear damping is present (h 2 = h = . . . = h n = 0) . In this way we can also 
compare an approximate formula that is given by Murphy in reference 7 , which, 
in our notation, is 


H Oe _ m Q + (3A)mg 

H ° m Q + (9/8)1112 


(85) 


For the sake of comparison, we write the corresponding case of (84) with 
fractions as 

Hoe = Kip + (5/8)1112 . . 

H o m Q + ( 15/l6 )m 2 

Both of these formulas have the correct limit as m 2 00 (Ho e /H 0 = 2/3) ) 
but when m Q = 1, formula (85) does not quite have the correct behavior as 
m 2 “*■ 0, yielding Hq /H q = 1 - (3/8)m2 + 0(m 2 2 ) instead of the exact value 
H 0e /H 0 = 1 - (5/l6)m 2 + 0(m 2 2 ) of the present analysis. For the stable-stable 
moment the approximate formula (86) is practically indistinguishable from the 
exact solution when plotted on the same graph, and (85) is also a good approx- 
imation for this case. These results are shown in figure 6(a). The main 
deviations occur in the stable-unstable and unstable -stable moments near the 
singular points, which we illustrate in figures 6(b) and (c). The approximate 
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The next check that we can make on the accuracy of equation (84) is to 
allow h 2 to be nonzero. A number of different cases were investigated, and 
one example for each type of linear-cubic moment is shown in figure 7* 

Figure 7(a) shows results for a stable- stable moment. The magnitude of the 
nonlinearities being considered on this plot can be visualized by observing 



(a) Stable- stable 1-3 static moment; stable- 
unstable or unstable -stable 1-3 damping 
moment . 



(b) Stable -unstable 1-3 static moment; stable- 
unstable or unstable- stable 1-3 damping 
moment . 



(c) Unstable- stable 1-3 static moment; stable - 
unstable or unstable- stable 1-3 damping 
moment . 


Figure 7*- Comparison of present analysis with 
numerical integrations . 






■'the small-scale plots of static and damping moments that appear. These are 
accurate representations of the nonlinearities, not sketches, and indicate a 
•wide departure from a linear system. The first thing to note in this figure is 
that equation ( 84 ) agrees to this scale exactly with equation (7l) or ( 74 ), 
which in turn agrees exactly with numerical integration of equation (l). The 
second thing to note is that a negative h 2 /m 2 (recall that the parameter 
h 2 /m 2 is a constant, independent of amplitude , for a given oscillator) can he 
obtained in two ways, H 0 negative and H 2 positive or Hp positive and H 2 
negative. The first case leads to a limit cycle, the second at a high enough 
amplitude to an unstable situation. The surprising thing is that Hp e /H 0 is 
the same for both of these cases, a fact that was revealed in the analysis but 
is not intuitively obvious. 

Also indicated in figure 7 (a) is the asymptotic behavior of H 0e /H 0 as 
given by equation ( 75 ) • It can be shown that expanding equation ( 84 ) for large 
m 2 (with h 4 ,h 6 , . . . = 0) yields the same values as equation (75) > except 
that the coefficient of mQh 2 /m 2 is O.289 instead of the more exact value 

O.296. 

Figures 7 (b) and (c) show representative cases for the stable-unstable and 
unstable- stable static moments. Once again the magnitude of the nonlinearities 
being covered is shown by the small-scale plots. For both of these cases, the 
accuracy of equation ( 84 ) is not as good as it was for the stable-stable 
moment, but the agreement is still excellent except in the immediate vicinity 
of the singular points. The final point to be made is that the curves for the 
limit amplitude, (-h 2 )j, presented in figure 5, can be closely duplicated by 
equation ( 84 ) except, again, in the immediate vicinity of the singular points. 


Comparison for a 1-9 Moment 

A more extreme nonlinear case was next investigated to check the adequacy 
of the approximate solution (82) . A 1-9 static moment and a 1-3 damping 
moment (MqjM^Hq,^ ^ 0 ) were considered. For this case, equation (82) 
appears as 


H o e = °o + O.30l6m e + m Q + 0.3596m e 
^o m Q + 0 . 904 7 m g m Q + 0 . 8940 m 8 

In this case, the closed-form solution (69) involves hyperelliptic integrals 
and cannot be evaluated by simple means, but a comparison can be made with 
results of numerical integrations of equation (l) . These comparisons are made 
in figure 8 for a stable- stable and unstable- stable static moment. (The 
stable-unstable case is not too interesting in that the moment is essentially 
linear almost to the singular point and then suddenly goes unstable.) Instead 
of the parameter h 2 /m 2 , which was pertinent for the cubic moment, the param- 
eter that remains constant for a given oscillator is now h 2 3 1^2 1 / m 8 * I6 e 
comparisons show that the approximate solution is not so good as it was for a 
linear-cubic static moment. However, for the stable-stable case, agreement is 
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still excellent. For the unstable -stable case, agreement is satisfactory 
except close to the singular point (i.e., m 8 = 5) • 




(a) Stable -stable 1-9 static moment; 1-3 (b) Unstable- stable 1-9 static moment; 1-3 

damping moment. damping moment. 

Figure 8.- Comparison of present analysis with 
numerical integrations. 


Comparison for a 1-3-5 Moment 

The final comparison that was made was for a l-3 _ 5 static moment in con- 
junction with several different damping moments, the most nonlinear of which 
was also 1-3-5* For this case, equation (82) appears as 

H° e _ m o + 0.6250m 2 + 0.4583m 4 mo + 0.6797m 2 + 0.5213m 4 

Ho m Q + 0.9375m 2 + 0 .9l67m 4 mo + 0.9297m 2 + 0.9067m 4 


0.7500m o + 0.5250m 2 + 0.4ll8m 4 
m Q + 0.9l88m 2 + 0.8929m 4 


( 88 ) 


Too many possibilities exist to cover this case in general. Hence, a particu- 
lar 1-3-5 static moment was selected as being very nonlinear, a stable- stable- 
unstable moment. This moment is shown in figure 9(a)* Three different 
damping moments were then chosen, consisting of the following: (l) linear 

term only, (2) linear-cubic, and (3) linear-cubic-quintic . The three damping 
moments are shown in figure 9(b)* Note that they were selected so that the 
average value of the damping moment between 0° and 20° was the same in each 
case . 
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.06 





y» d eg 

(a) 1-3-5 static moment. 



(h) Damping moments. 


Figure 9*- Moments used for more complex 
comparisons. 



y 0 ,deg 


Figure 10.- Comparison of present analysis with 
numerical integrations; 1-3-5 static moment; 
various damping moments. 


Equation (88) is compared with 
numerical 'integrations of equation (l) 
for the three cases in figure 10. To 
the scale plotted and for the range 
covered, no differences existed. 

All of these comparisons that 
have been presented indicate that 
equation (82) is surprisingly accurate 
over a wide range of nonlinearities in 
both the static and damping moment. 
With this established, we will next 
concern ourselves with determination 
of the coefficients M s and H s from 
data. 


^3 




determination of nonlinear parameters from data 


Formula (82) with (83) can be useful in determining the nonlinear damping 
parameters of an oscillator from a set of observed oscillations. Let us assume 
that a given oscillator is governed by the differential equation (l) . We wish 
to infer from observed oscillations the appropriate values of H Q ,H 2 ,H 4 , . . . 
and M Q ,M 2 ,M 4 , .... 

We assume that the damping is small so that h = H Q /v « 1 (or more 
appropriately, since v = |Mny a n | 1/2 is actually arbitrary, we could also say 
H 0 /w « 1 , where w is the frequency and is proportional to v) . This assump- 
tion is normally satisfied in ballistic-range testing. In addition, we assume 
that the frequency and maximum amplitudes for each half cycle of the motion 
can be accurately determined from the data* The effective linear damping for 
each half cycle is then computed by 


H o e s I “ 


r n 


r n+i 


(89) 


Using an approximation for the frequency quadrature (eq. (67b) ) that is 
developed in reference 13* we can write an expression for the frequency as a 
function of the parameters ... as 


where 


GO 


= M q + 


^ A n M ny a n 


n>o 


(90a) 



(90b) 


This formula is a good approximation for large and small nonlinearities, except 
in the vicinity of the singular points of static instability. The understand- 
ing and use_of equation (90a) makes it possible to extract the appropriate 
values of M n by fitting equation ( 90 a) to a set of frequency versus amplitude 
data that has been measured for a given oscillator. For instance, if an 
oscillator is governed by a linear-cubic static moment, then equation (90) 
yields 

w 2 = M 0 + ~ M 2 y a 2 


The frequency squared then gives a straight line when plotted against y a 2 . 

The slope of the curve yields the_value of M 2 and the w 2 intercept yields 
M q . When h « 1, we note that M n = M^l + 0 (h 2 )]. Hence, for small damping, 
the values of the actual static parameters (M n ) are determined. The cubic 
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moment is the simplest case; for more general cases the static moment polyno- 
mial (eq. ( 90 a)) is terminated so that an adequate fit of the frequency versus 
maximum amplitude data is obtained (ref. 13 ) - 


A similar but slightly more complicated situation exists for determining 
the damping parameters H n . In this case, the static parameters play an impor- 
tant role, whereas for the determination of the static-moment parameters Mq, 
the effect of the damping parameters can be neglected. The formula for Hq 6 
as a function of • • . can be written as 


H 0e = Bo H o + £ B nHny a n 
n>o 

where equation ( 69 ) is used to determine B n as 



or, from equations (82) and ( 83 ) , approximately as 


1 s>o 


Y, ^s^a 5 


®n 2 + n 


8 M y s 
H ns s^a 


s>o 


( 91 a) 


(91b) 


(91c) 


Expression (91a) is analogous in form to (20a). In this case, however, the 
coefficients B n are functions of M Q ,M 1 ,M 2 , . . . and y a , whereas the coef- 
ficients A n are constants. In the determination of the damping, however, the 
static parameters can be regarded as already determined by the use of ( 90 ) . 
Hence, the coefficients B n are functions of y a insofar as determination of 
the damping is concerned. 

It is now useful to divide both sides of equation (91a) by B Q , obtaining 




n>o 


(92) 
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Because B 0 >B 1 ,B 2 , • • • each vary with y a in a similar manner, the ratios 
B n /B 0 vary slowly with y a . Equation (92) can now be used to study a set of * 
damping data in a manner analogous to studying frequency data with 
equation (90) . 


As a simple example, assume that a set of data for a given oscillator has 
been obtained and that the static moment, completely arbitrary, has been deter- 
mined with the use of formula (90) . Let it be assumed that the damping moment 
can be described by a linear term and a single arbitrary nonlinear term so 
that (92) appears as 

^ - H o + H n y a n ( 93 ) 

■°o 

where n is not known. Choosing a particular value of n determines B n /B 0 . 
Then one can plot E Qq /b 0 versus (B n /B 0 )y a n . When the right value of n is 
chosen, the data will fall on a straight line. The slope of the line yields 
H n , and the intercept of the line with the H 0e /B 0 axis yields H Q . This will 
he true regardless of the particular form of the static moment. 

For an arbitrary damping moment described by a given set of data over an 
amplitude range, there may be several or many combinations of the parameters 
Ho^H x ,H 2 , . . . that will fit the data. This situation also occurs in deter- 
mining the static-moment parameters (ref. 13)- One must often settle for a 
member of a class of moments that gives a good fit. Thus, a certain amount of 
experience may be useful in analyzing the data. For a set of data including 

small angles of y a as well as large, 
one should plot the damping data as 
Ho e /B 0 versus (B 2 /B 0 )y a 2 (or as a first 
try against y a 2 since B 2 /B 0 varies 
slowly with y a ) . If the data fall on 
a straight line, the damping moment is 
a cubic as described by (93) with 
n = 2. If the data deviate from a 
straight line, which might be expected 
to occur at larger values of y a , one 
can estimate the magnitude and sign of 
the next higher order damping term 
needed to fit the data. In such a man- 
ner the appropriate form of damping 
polynomial tends to suggest itself. 

The more data that are available, 
the better the damping moment can be 
defined. Finally the data can be 
fitted by some curve-fitting technique 
for various combinations of damping 
parameters. The combinations that fit 
the data best are the appropriate ones 
to represent the damping moment. 

'a 

In figure 11 the results of fig- 

Figure li.- Form of data to reveal type of 10 are Pitted in the manner sug- 

nonlinearity in dairying . gested. Results for the linear-cubic 
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damping moment fall on a straight line as expected. The deviation of the 
* 1 - 3-5 damping data from a straight line indicates that the damp ing moment is 
more complicated than linear cubic; the fact that the results curve downward 
from a straight line indicates that the next higher order term required is a 
destabilizing term, which is, in fact, the case. 


CONCLUDING REMARKS 


The preceding investigation has considered the effects of nonlinearities 
on danced oscillations. It was found that this investigation could be con- 
ducted analytically in two ways. With large damp ing nonlinearities and small 
"effective" static-moment nonlinearities, an approximate solution could be 
established by means of an integral equation. On the other hand, when the 
damping is small but the static nonlinearities arbitrarily large (a common sit- 
uation in aerodynamics), an approximate solution could be obtained by means of 
an equivalent first-order differential equation. 

When the damping is small and the static nonlinearities large, it is con- 
venient to study the characteristics of the nonlinear oscillations by means of 
a parameter called the "effective linear damping." This parameter is readily 
obtainable from experiment and offers a means of deducing the nonlinear damp- 
ing characteristics of an oscillator from a set of observed oscillations. 
Comparisons with some exact solutions and with numerical solutions showed that 
the approximate formula derived herein for the effective linear damping is 
extremely accurate over a wide range of nonlinearities. The limit cycle of a 
nonlinear oscillation occurs when the "effective linear damping" is zero, and 
the limiting amplitude is easy to obtain for various combinations of the non- 
linear damping and static moments. The concept of an "effective linear 
damping" is thus very useful in analyzing nonlinear damping. 

The integral equation ( 36 ) and the first-order equation (39) are believed 
to be novel formulations for nonlinear oscillation problems. It is possible 
that further investigation of these equations would be fruitful. 


Ames Research Center 

National Aeronautics and Space Administration 
Moffett Field, Calif., June 27, 19^6 
124-07-02-11-21 
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APPENDIX A 




MODIFIED EQUATIONS OF KRYLOFF AND BOGOLIUBOFF 


Here we wish to modify the approximation methods of Kryloff and 
Bogoliuboff (ref. 4) for nonlinear oscillations so that the pertinent equations 
fit conveniently into the framework of the present analysis. We shall forego 
the Fourier series expansions of the usual Kryloff -Bogoliuboff analysis and 
find forms of the amplitude and frequency equations that are amenable to appro- 
priate transformations that in turn yield a basic exact solution as a special 
case. 


Following Kryloff and Bogoliuboff, consider a nonlinear equation of the 
following form : 




(Al) 


where f is an arbitrary function of x and dx/dt. If f = 0, then the solu- 
tion is given by 


x = a cos(vt + cp) 


(A2) 


where a and 9 are arbitrary constants. For the general equation (Al) , we 
assume equation (A2) to be a solution with a = a(t) and 9 = 9(t) suitably 
determined functions of t . 

Let us choose a(t) and 9(t) by following the method of variation of 
parameters for linear systems. The derivative of equation (A2) is 


— = ^ cos(vt +9) - a (v + sin(vt + 
dt dt \ dt J 


9 ) 


(A3) 


Since we have two functions to determine, a(t) and 9(t) , we have two conditions 
at our disposal. One, of course, will be the original differential equa- 
tion (Al) . The other we choose as 


— cos(vt + 9) 



sin(vt + 9) = 0 


(A4) 


so that equation (A 3 ) now reads simply 


— = -av sin(vt + 9) (A 5 ) 

dt 
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just as for the case f = 0 . Differentiating equation (A 5 ) and substituting 
into the original equation (Al) , we get 


sin(vt + 9 ) + a cos(vt + 9 ) = i f[a cos(vt + 9 ), -av sin(vt + 9 )] 

(A 6 ) 

Expressions (A4) and (A 6 ) are two equations for da/dt and d 9 /dt. Solving 
for da/dt and d 9 /dt, we obtain 

i f[a cos(vt + 9 ), -av sin(vt + 9)]sin(vt + 9 ) (ATa) 

4J = — f[a cos(vt + 9 ), -av sin(vt + 9 )]cos(vt + 9 ) (A7b) 


These are the basic equations of Kryloff and Bogoliuboff . 

Equations (A7a, b) are coupled and thus must be solved simultaneously. 
Equation (A7a) can be uncoupled, however, if we introduce, in place of t, a 
new independent variable t defined as follows: 

T = vt + 9 (A 8 a) 


It follows that 


It 

dt 



(A8b) 


Equations (A7a, b) may now be written in terms of t as follows: 


da _ v 
dr 


f(a cos t, -av sin -r)sin t 


1 + — f(a cos t, -av sin t)cos t 
av^ 


(A9a) 


— f(a cos t, -av sin t)cos t 
d9 _ av^ 

dT 


1 + — p f(a cos t, -av sin t)cos t 


av 


(A9b) 


dT _ 
dt 


- V! 


1 + f(a cos t, -av sin t)cos t 


av 


(A9c) 


We shall term equations (A9a, b, c) the modified Kryloff -Bogoliuboff equa- 
tions. Equations (A9a) and (A9c) are the pertinent ones since the solution to 
equation (Al) is given by y = a(-r)eos t* Equation (A9a) is a first-order 
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